\section{Conclusion}
Most existing methods in LiDAR-based 3D detection rely on hand-crafted feature representations, for example, a bird's eye view projection. In this paper, we remove the bottleneck of manual feature engineering and propose VoxelNet, a novel end-to-end trainable deep architecture for point cloud based 3D detection. Our approach can operate directly on sparse 3D points and capture 3D shape information effectively. We also present an efficient implementation of VoxelNet that benefits from point cloud sparsity and parallel processing on a voxel grid. Our experiments on the KITTI car detection task show that VoxelNet outperforms state-of-the-art LiDAR based 3D detection methods by a large margin. On more challenging tasks, such as 3D detection of pedestrians and cyclists, VoxelNet also demonstrates encouraging results showing that it provides a better 3D representation. Future work includes extending VoxelNet for joint LiDAR and image based end-to-end 3D detection to further improve detection and localization accuracy.

\section{Training Details}
\label{sec:training_details}
In this section, we explain the implementation details of the VoxelNet and the training procedure.
\subsection{Network Details}
Our experimental setup is based on the LiDAR specifications of the KITTI dataset~\cite{REF:Geiger2012CVPR}.

\noindent{\bf Car Detection }  For this task, we consider point clouds within the range of $[-3, 1] \times [-40, 40] \times [0, 70.4]$ meters along Z, Y, X axis respectively. Points that are projected outside of image boundaries are removed~\cite{REF:cvpr17chen}. We choose a voxel size of $v_D=0.4, v_H=0.2, v_W=0.2$ meters, which leads to $D'=10$, $H'=400$, $W'=352$. We set $T=35$ as the maximum number of randomly sampled points in each non-empty voxel. We use two VFE layers $\textrm{VFE-1}(7, 32)$ and $\textrm{VFE-2}(32, 128)$. The final FCN maps VFE-2 output to $\mathbb{R}^{128}$. Thus our feature learning net generates a sparse tensor of shape $128 \times 10 \times 400 \times 352$. To aggregate voxel-wise features, we employ three convolution middle layers sequentially as Conv3D(128, 64, 3, (2,1,1), (1,1,1)), Conv3D(64, 64, 3, (1,1,1), (0,1,1)), and Conv3D(64, 64, 3, (2,1,1), (1,1,1)), which yields a 4D tensor of size $64 \times 2 \times 400 \times 352$. After reshaping, the input to RPN is a feature map of size $128 \times 400 \times 352$, where the dimensions correspond to channel, height, and width of the 3D tensor. Figure~\ref{fig:RPN} illustrates the detailed network architecture for this task. Unlike \cite{REF:cvpr17chen}, we use only one anchor size, $l^a=3.9, w^a=1.6, h^a=1.56$ meters, centered at $z_c^a = -1.0$ meters with two rotations, 0 and 90 degrees.
Our anchor matching criteria is as follows: An anchor is considered as positive if it has the highest Intersection over Union (IoU) with a ground truth or its IoU with ground truth is above 0.6 (in bird's eye view). An anchor is considered as negative if the IoU between it and all ground truth boxes is less than 0.45. We treat anchors as don't care if they have $0.45 \leq \textrm{IoU} \leq 0.6$ with any ground truth. We set $\alpha = 1.5$ and $\beta = 1$ in Eqn.~\ref{eqn:loss_function}.

\noindent{\bf Pedestrian and Cyclist Detection } The input range\footnote{Our empirical observation suggests that beyond this range, LiDAR returns from pedestrians and cyclists become very sparse and therefore detection results will be unreliable.} is $[-3, 1] \times [-20, 20] \times [0, 48]$ meters along Z, Y, X axis respectively. We use the same voxel size as for car detection, which yields $D=10$, $H=200$, $W=240$. We set $T=45$ in order to obtain more LiDAR points for better capturing shape information. The feature learning network and convolutional middle layers are identical to the networks used in the car detection task. For the RPN, we make one modification to block 1 in Figure~\ref{fig:RPN} by changing the stride size in the first 2D convolution from 2 to 1. This allows finer resolution in anchor matching, which is necessary for detecting pedestrians and cyclists. We use anchor size $l^a=0.8, w^a=0.6, h^a=1.73$ meters  centered at $z_c^a = -0.6$ meters with 0 and 90 degrees rotation for pedestrian detection and use anchor size $l^a=1.76, w^a=0.6, h^a=1.73$ meters centered at $z_c^a = -0.6$ with 0 and 90 degrees rotation for cyclist detection.
The specific anchor matching criteria is as follows: We assign an anchor as postive if it has the highest IoU with a ground truth, or its IoU with ground truth is above 0.5. An anchor is considered as negative if its IoU with every ground truth is less than 0.35. For anchors having $0.35 \leq \textrm{IoU} \leq 0.5$ with any ground truth, we treat them as don't care.

During training, we use stochastic gradient descent (SGD) with learning rate 0.01 for the first 150 epochs and decrease the learning rate to 0.001 for the last 10 epochs. We use a batchsize of 16 point clouds.

\subsection{Data Augmentation}
With less than 4000 training point clouds, training our network from scratch will inevitably suffer from overfitting. To reduce this issue, we introduce three different forms of data augmentation. The augmented training data are generated on-the-fly without the need to be stored on disk~\cite{REF:NIPS2012_4824}.

Define set $\mathbf{M} = \{ \mathbf{p}_i=[x_i, y_i, z_i, r_i]^T \in \mathbb{R}^4 \}_{i=1, \ldots, N}$ as the whole point cloud, consisting of $N$ points.
We parameterize a  3D bouding box $\mathbf{b}_i$ as $(x_c, y_c, z_c, l, w, h, \theta)$, where $x_c, y_c, z_c$ are center locations, $l, w, h$ are length, width, height, and $\theta$ is the yaw rotation around Z-axis. We define $\Omega_i=\{\mathbf{p} | x \in [x_c - l/2, x_c + l/2], y \in [y_c - w/2, y_c + w/2], z \in [z_c - h/2, z_c + h/2], \mathbf{p} \in \mathbf{M}\}$ as the set containing all LiDAR points within $\mathbf{b}_i$, where $\mathbf{p}=[x, y, z, r]$ denotes a particular LiDAR point in the whole set $\mathbf{M}$.

The first form of data augmentation applies perturbation independently to each ground truth 3D bounding box together with those LiDAR points within the box. Specifically, around Z-axis we rotate $\mathbf{b}_i$ and the associated $\Omega_i$ with respect to $(x_c, y_c, z_c)$ by a uniformally distributed random variable $\Delta\theta \in [-\pi/10, +\pi/10]$. Then we add a translation $(\Delta x, \Delta y, \Delta z)$ to the XYZ components of $\mathbf{b}_i$ and to each point in $\Omega_i$, where $\Delta x$, $\Delta y$, $\Delta z$ are drawn independently from a Gaussian distribution with mean zero and standard deviation 1.0. To avoid physically impossible outcomes, we perform a collision test between any two boxes after the perturbation and revert to the original if a collision is detected. Since the perturbation is applied to each ground truth box and the associated LiDAR points independently, the network is able to learn from substantially more variations than from the original training data.

Secondly, we apply global scaling to all ground truth boxes $\mathbf{b}_i$ and to the whole point cloud $\mathbf{M}$. Specifically, we multiply the XYZ coordinates and the three dimensions of each $\mathbf{b}_i$, and the XYZ coordinates of all points in $\mathbf{M}$ with a random variable drawn from uniform distribution $[0.95, 1.05]$. Introducing global scale augmentation improves robustness of the network for detecting objects with various sizes and distances as shown in image-based classification~\cite{REF:Simonyan14c, REF:journals/corr/Howard13} and detection tasks~\cite{REF:Girshick:2015, REF:ResNet2016}.

Finally, we apply global rotation to all ground truth boxes $\mathbf{b}_i$ and to the whole point cloud $\mathbf{M}$. The rotation is applied along Z-axis and around $(0, 0, 0)$. The global rotation offset is determined by sampling from uniform distribution $[-\pi/4, +\pi/4]$. By rotating the entire point cloud, we simulate the vehicle making a turn.
\section{Discussion}
\section{Experiments}
\label{sec:Experiment}
We evaluate VoxelNet on the KITTI 3D object detection benchmark~\cite{REF:Geiger2012CVPR} which contains 7,481 training images/point clouds and 7,518 test images/point clouds, covering three categories: \textit{Car}, \textit{Pedestrian}, and \textit{Cyclist}. For each class, detection outcomes are evaluated based on three difficulty levels: \textit{easy}, \textit{moderate}, and \textit{hard}, which are determined according to the object size, occlusion state, and truncation level. Since the ground truth for the test set is not available and the access to the test server is limited, we conduct comprehensive evaluation using the protocol described in~\cite{REF:nips15chen,REF:cvpr16chen,REF:cvpr17chen} and subdivide the training data into a training set and a validation set, which results in 3,712 data samples for training and 3,769 data samples for validation. The split avoids samples from the same sequence being included in both the training and the validation set~\cite{REF:cvpr16chen}. Finally we also present the test results using the KITTI server.

For the \textit{Car} category, we compare the proposed method with several top-performing algorithms, including image based approaches: Mono3D~\cite{REF:cvpr16chen} and 3DOP~\cite{REF:nips15chen}; LiDAR based approaches: VeloFCN~\cite{REF:VeloFCN} and 3D-FCN~\cite{REF:3DFCN}; and a multi-modal approach MV~\cite{REF:cvpr17chen}. Mono3D~\cite{REF:cvpr16chen}, 3DOP~\cite{REF:nips15chen} and MV~\cite{REF:cvpr17chen}  use a pre-trained model for initialization whereas we train VoxelNet from scratch using only the LiDAR data provided in KITTI.

To analyze the importance of end-to-end learning, we implement a strong baseline that is derived from the VoxelNet architecture but uses hand-crafted features instead of the proposed feature learning network. We call this model the hand-crafted baseline (HC-baseline). HC-baseline uses the bird's eye view features described in
\cite{REF:cvpr17chen} which are computed at $0.1$m resolution. Different from \cite{REF:cvpr17chen}, we increase the number of height channels from 4 to 16 to capture more detailed shape information-- further increasing the number of height channels did not lead to performance improvement. We replace the convolutional middle layers of VoxelNet with similar size 2D convolutional layers, which are Conv2D(16, 32, 3, 1, 1), Conv2D(32, 64, 3, 2, 1), Conv2D(64, 128, 3, 1, 1). Finally RPN is identical in VoxelNet and HC-baseline. The total number of parameters in HC-baseline and VoxelNet are very similar. We train the HC-baseline using the same training procedure and data augmentation described in Section~\ref{sec:training_details}.

\begin{figure*}[!t]
\centering
    \begin{tabular}{ccc}
    \centering

\vspace{-0.30cm}

        \hspace{-0.75cm} \includegraphics[width=0.38\linewidth]{001755_vox_rgb.png} & \hspace{-1.25cm}

        \includegraphics[width=0.38\linewidth]{005467_ped.png} & \hspace{-1.25cm}

        \includegraphics[width=0.38\linewidth]{006407_cyclist.png} \\

\vspace{-0.30cm}

        \hspace{-0.75cm} \includegraphics[width=0.38\linewidth]{002254_vox_rgb.png} & \hspace{-1.25cm}

        \includegraphics[width=0.38\linewidth]{005840_ped.png} & \hspace{-1.25cm}

        \includegraphics[width=0.38\linewidth]{003763_cyclist.png} \\

\vspace{-0.30cm}

        \hspace{-0.75cm} \includegraphics[width=0.38\linewidth]{000519_car.png} & \hspace{-1.25cm}

        \includegraphics[width=0.38\linewidth]{006210_ped.png} & \hspace{-1.25cm}

        \includegraphics[width=0.38\linewidth]{000383_cyclist.png} \\

\vspace{-0.30cm}

        \hspace{-0.75cm} \includegraphics[width=0.38\linewidth]{000669_car.png} & \hspace{-1.25cm}

        \includegraphics[width=0.38\linewidth]{005910_ped.png} & \hspace{-1.25cm}

        \includegraphics[width=0.38\linewidth]{003677_cyclist.png} \\

        \hspace{-0.75cm} Car & \hspace{-1.25cm} Pedestrian & \hspace{-1.25cm} Cyclist \\
    \end{tabular}

    \caption{Qualitative results. For better visualization  3D boxes detected using LiDAR are projected on to the RGB images.}
    \label{fig:Qualitative_demo}
\end{figure*}

\subsection{Evaluation on KITTI Validation Set}
\paragraph{Metrics}
We follow the official KITTI evaluation protocol, where the IoU threshold is 0.7 for class \textit{Car} and is 0.5 for class \textit{Pedestrian} and \textit{Cyclist}. The IoU threshold is the same for both bird's eye view and full 3D evaluation. We compare the methods using the average precision (AP) metric.

\paragraph{Evaluation in Bird's Eye View}

The evaluation result is presented in Table~\ref{table:BEV_detection_comparison}. VoxelNet consistently outperforms all the competing approaches across all three difficulty levels. HC-baseline also achieves satisfactory performance compared to the state-of-the-art~\cite{REF:cvpr17chen}, which shows that our base region proposal network (RPN) is effective. For \textit{Pedestrian} and \textit{Cyclist} detection tasks in bird's eye view, we compare the proposed VoxelNet with HC-baseline. VoxelNet yields substantially higher AP than the HC-baseline for these more challenging categories, which shows that end-to-end learning is essential for point-cloud based detection.

We would like to note that \cite{REF:3DFCN} reported 88.9\%, 77.3\%, and 72.7\% for \textit{easy}, \textit{moderate}, and \textit{hard} levels respectively, but these results are obtained based on a different split of 6,000 training frames and $\sim$1,500 validation frames, and they are not directly comparable with algorithms in Table~\ref{table:BEV_detection_comparison}. Therefore, we do not include these results in the table.

\paragraph{Evaluation in 3D}
Compared to the bird's eye view detection, which requires only accurate localization of objects in the 2D plane, 3D detection is a more challenging task as it requires finer localization of shapes in 3D space. Table~\ref{table:3D_detection_comparison} summarizes the  comparison. For the class \textit{Car}, VoxelNet significantly outperforms all other approaches in AP across all difficulty levels. Specifically, using only LiDAR, VoxelNet significantly outperforms the state-of-the-art method MV (BV+FV+RGB)~\cite{REF:cvpr17chen} based on LiDAR+RGB, by 10.68\%, 2.78\% and 6.29\% in \textit{easy}, \textit{moderate}, and \textit{hard} levels respectively. HC-baseline achieves similar accuracy to the MV~\cite{REF:cvpr17chen} method.

As in the bird's eye view evaluation, we also compare VoxelNet with HC-baseline on 3D \textit{Pedestrian} and \textit{Cyclist} detection. Due to the high variation in 3D poses and shapes, successful detection of these two categories requires better 3D shape representation. As shown in Table~\ref{table:3D_detection_comparison} the improved performance of VoxelNet is emphasized for more challenging 3D detection tasks (from $\sim$8\% improvement in bird's eye view to $\sim$12\% improvement on 3D detection) which suggests that VoxelNet is more effective in capturing 3D shape information than hand-crafted features.

\subsection{Evaluation on KITTI Test Set}

We evaluated VoxelNet on the KITTI test set by submitting detection results to the official server. The results are summarized in Table~\ref{table:test_set_result}. VoxelNet,  significantly outperforms the previously published state-of-the-art~\cite{REF:cvpr17chen}  in all the tasks (bird's eye view and 3D detection) and all difficulties. We would like to note that many of the other leading methods listed in KITTI benchmark  use both RGB images and LiDAR point clouds whereas VoxelNet  uses only LiDAR.

We present several 3D detection examples in Figure~\ref{fig:Qualitative_demo}. For better visualization  3D boxes detected using LiDAR are projected on to the RGB images. As shown, VoxelNet provides highly accurate 3D bounding boxes in all categories.

The inference time for the VoxelNet is 225ms where the voxel input feature computation takes 5ms, feature learning net takes 20ms, convolutional middle layers take 170ms, and region proposal net takes 30ms on a TitanX GPU and 1.7Ghz CPU.

\begin{table}[!t]
\centering
\scalebox{0.9}{
\begin{tabular}{| c | c | c | c |}
\hline
{\bf Benchmark} & {\bf Easy} & {\bf Moderate} & {\bf Hard} \\
\hline
Car (3D Detection) & 77.47  & 65.11  & 57.73 \\
Car (Bird's Eye View) & 89.35  & 79.26  & 77.39 \\
\hline
Pedestrian (3D Detection) & 39.48  & 33.69  & 31.51 \\
Pedestrian (Bird's Eye View) & 46.13  & 40.74  & 38.11 \\
\hline
Cyclist (3D Detection) & 61.22  & 48.36  & 44.37 \\
Cyclist (Bird's Eye View) & 66.70  & 54.76  & 50.55  \\
\hline
\end{tabular}
}
\caption{Performance evaluation on KITTI test set.}
\label{table:test_set_result}
\end{table}

\section{Introduction}

\begin{figure}[!ht]
\centering
    \includegraphics[width=0.9\linewidth]{1st_figure_cartoon_new.png}
\caption{VoxelNet directly operates on the raw point cloud (no need for feature engineering)  and produces the 3D detection results using a single end-to-end trainable network. }
\label{fig:1st_figure}
\end{figure}

\begin{figure*}[!th]
\centering
    \includegraphics[width=0.80\linewidth]{VoxelNet_new.png}
    \vspace{-0.1cm}
\caption{VoxelNet architecture. The feature learning network takes a raw point cloud as input, partitions the space into voxels, and transforms points within each voxel to a vector representation characterizing the shape information. The space is represented as a sparse 4D tensor. The convolutional middle layers processes the 4D tensor to aggregate spatial context. Finally, a RPN generates the 3D detection.}
\label{fig:VoxelNet_flowchart}
\end{figure*}

Point cloud based 3D object detection is an important component of a variety of real-world applications, such as autonomous navigation~\cite{REF:Geiger2012CVPR,REF:Drone}, housekeeping robots~\cite{REF:Housekeeping}, and augmented/virtual reality~\cite{REF:AR}. Compared to image-based detection, LiDAR provides reliable depth information that can be used to accurately localize objects and characterize their shapes~\cite{REF:3DFCN, REF:cvpr17chen}. However, unlike images,  LiDAR point clouds are sparse and have highly variable point density, due to factors such as non-uniform sampling of the 3D space, effective range of the sensors, occlusion, and the relative pose. To handle these challenges, many approaches manually crafted feature representations for point clouds that are tuned for 3D object detection. Several methods project point clouds into a perspective view and apply image-based feature extraction techniques~\cite{REF:FusionDPM-IROS14, REF:MV-RGBD-RF2015,REF:VeloFCN}. Other approaches rasterize point clouds into a 3D voxel grid and encode each voxel with hand-crafted features~\cite{REF:Wang-RSS-15,REF:Vote3Deep,REF:Song2014,REF:DeepSlidingShapes,REF:3DFCN,REF:cvpr17chen}. However, these manual design choices introduce an information bottleneck that prevents these approaches from effectively exploiting 3D shape information and the required invariances for the detection task.  A major breakthrough in recognition~\cite{REF:NIPS2012_4824} and detection~\cite{REF:girshick2014rich} tasks on images was due to moving from hand-crafted features to machine-learned features.

Recently, Qi \textit{et al.}\cite{REF:qi2017pointnet} proposed PointNet, an end-to-end deep neural network that learns point-wise features directly from point clouds. This approach demonstrated impressive results on 3D object recognition, 3D object part segmentation, and point-wise semantic segmentation tasks. In \cite{REF:qi2017pointnetplusplus}, an improved version of PointNet was introduced which enabled the network to learn local structures at different scales. To achieve satisfactory results, these two approaches trained feature transformer networks on all input points ($\sim$1k points). Since typical point clouds obtained using LiDARs contain $\sim$100k points, training the architectures as in~\cite{REF:qi2017pointnet,REF:qi2017pointnetplusplus} results in high computational and memory requirements.   Scaling up 3D feature learning networks to orders of magnitude more points and to 3D detection tasks are the main challenges that we address in this paper.

Region proposal network (RPN)~\cite{REF:NIPS2015_5638} is a highly optimized algorithm for efficient object detection~\cite{REF:ResNet2016,REF:cvpr17chen,REF:YOLO9000_2017,REF:SSD_Liu2016}. However, this approach requires data to be dense and organized in a tensor structure (e.g. image, video) which is not the case for typical LiDAR point clouds. In this paper, we close the gap between point set feature learning and RPN for 3D detection task.

We present VoxelNet, a generic 3D detection framework that simultaneously learns a discriminative feature representation from point clouds and predicts accurate 3D bounding boxes, in an end-to-end fashion, as shown in Figure~\ref{fig:VoxelNet_flowchart}.
We design a novel voxel feature encoding (VFE) layer, which enables inter-point interaction within a voxel, by combining point-wise features with a locally aggregated feature. Stacking multiple VFE layers allows learning complex features for characterizing local 3D shape information. Specifically, VoxelNet divides the point cloud into equally spaced 3D voxels, encodes each voxel via stacked VFE layers, and then 3D convolution further aggregates local voxel features, transforming the point cloud into a high-dimensional volumetric representation. Finally, a RPN consumes the volumetric representation and yields the detection result. This efficient algorithm benefits both from the sparse point structure and efficient parallel processing on the voxel grid.

We evaluate VoxelNet on the bird's eye view detection and the full 3D detection tasks, provided by the KITTI benchmark~\cite{REF:Geiger2012CVPR}. Experimental results show that VoxelNet outperforms  the state-of-the-art LiDAR based 3D detection methods by a large margin. We also demonstrate that VoxelNet achieves highly encouraging results in detecting pedestrians and cyclists from LiDAR point cloud.

\subsection{Related Work}
\label{sec:RelatedWork}

Rapid development of 3D sensor technology has motivated researchers to develop efficient representations to detect and localize objects in point clouds. Some of the earlier methods for feature representation are \cite{REF:StructuralIndexing_Medioni92,REF:COSMOS_Jain97,REF:PointSignatures_Chua1997,REF:SpinImageHebert99,REF:Tuzel2014,REF:FPFH_Beetz2009,REF:RoboticAssembly_Oncel2012,REF:Mian2010,REF:Nishino2010,REF:MS_POSE_Blake2011,REF:bo_iros11}. These hand-crafted features yield satisfactory results when rich and detailed 3D shape information is available. However their inability to adapt to more complex shapes and scenes, and learn required invariances from data resulted in limited success for uncontrolled scenarios such as autonomous navigation.

Given that images provide detailed texture information, many algorithms infered the 3D bounding boxes from 2D images~\cite{REF:nips15chen,REF:cvpr16chen,REF:xiang_cvpr15,REF:Zia2013,REF:Zia2014,REF:SFM2015}. However, the accuracy of  image-based 3D detection approaches are bounded by the accuracy of the depth estimation.

Several LIDAR based 3D object detection techniques utilize a voxel grid representation.  \cite{REF:Wang-RSS-15, REF:Vote3Deep} encode each nonempty voxel with 6 statistical quantities that are derived from all the points contained within the voxel. \cite{REF:Song2014} fuses multiple local statistics to represent each voxel. \cite{REF:DeepSlidingShapes} computes the truncated signed distance on the voxel grid. \cite{REF:3DFCN} uses binary encoding for the 3D voxel grid. \cite{REF:cvpr17chen} introduces a multi-view representation for a LiDAR point cloud by computing a multi-channel feature map in the bird's eye view and the cylindral coordinates in the frontal view. Several other studies project point clouds onto a perspective view and then use image-based feature encoding schemes~\cite{REF:FusionDPM-IROS14, REF:MV-RGBD-RF2015,REF:VeloFCN}.

There are also several multi-modal fusion methods that combine images and LiDAR to improve detection accuracy~\cite{REF:Enzweiler2011,REF:Gonzalez2017,REF:cvpr17chen}. These methods provide improved performance compared to LiDAR-only 3D detection, particularly for small objects (pedestrians, cyclists) or when the objects are far, since cameras provide an order of magnitude more measurements than LiDAR. However the need for an additional camera that is time synchronized and calibrated with the LiDAR restricts their use and makes the solution more sensitive to sensor failure modes. In this work we focus on LiDAR-only detection.

\subsection{Contributions}

\begin{itemize}
	\item  We propose a novel end-to-end trainable deep architecture for point-cloud-based 3D detection, VoxelNet, that directly operates on sparse 3D points and avoids information bottlenecks introduced by manual feature engineering.
	\item We present an efficient method to implement VoxelNet which benefits both from the sparse point structure and efficient parallel processing on the voxel grid.
	\item We conduct experiments on KITTI benchmark and show that VoxelNet produces state-of-the-art results in LiDAR-based car, pedestrian, and cyclist detection benchmarks.
\end{itemize}

\section{VoxelNet}
\label{sec:VoxelNet}
In this section we explain the architecture of VoxelNet, the loss function used for training, and an efficient algorithm to implement the network.

\subsection{VoxelNet Architecture}
The proposed VoxelNet consists of three functional blocks: (1) Feature learning network, (2) Convolutional middle layers, and (3) Region proposal network~\cite{REF:NIPS2015_5638}, as illustrated in Figure~\ref{fig:VoxelNet_flowchart}. We provide a detailed introduction of VoxelNet in the following sections.

\subsubsection{Feature Learning Network}
\label{subsec:Feature_Learning}
\noindent{\bf Voxel Partition } Given a point cloud, we subdivide the 3D space into equally spaced voxels as shown in Figure~\ref{fig:VoxelNet_flowchart}. Suppose the point cloud encompasses 3D space with range $D$, $H$, $W$ along the Z, Y, X axes respectively. We define each voxel of size $v_D$, $v_H$, and $v_W$ accordingly. The resulting 3D voxel grid is of size $D' = D/v_D, H' = H/v_H, W' = W/v_W$. Here, for simplicity, we assume $D$, $H$, $W$ are a multiple of  $v_D$, $v_H$, $v_W$.

\noindent{\bf Grouping } We group the points according to the voxel they reside in. Due to factors such as distance, occlusion, object's relative pose, and non-uniform sampling, the LiDAR point cloud is sparse and has highly variable point density throughout the space. Therefore, after grouping, a voxel will contain a variable number of points. An illustration is shown in Figure~\ref{fig:VoxelNet_flowchart}, where Voxel-1 has significantly more points than Voxel-2 and Voxel-4, while Voxel-3  contains no point.

\noindent{\bf Random Sampling } Typically a high-definition LiDAR point cloud is composed of $\sim$100k points. Directly processing all the points not only imposes increased memory/efficiency burdens on the computing platform, but also highly variable point density throughout the space might bias the detection. To this end, we  randomly sample a fixed number, $T$, of points from those voxels containing more than $T$ points. This sampling strategy has two purposes, (1) computational savings (see Section~\ref{subsec:implementation} for details); and (2) decreases the imbalance of points between the voxels which reduces the sampling bias, and adds more variation to training.

\noindent{\bf Stacked Voxel Feature Encoding } The key innovation is the chain of VFE  layers. For simplicity, Figure~\ref{fig:VoxelNet_flowchart} illustrates the hierarchical feature encoding process for one voxel. Without loss of generality, we use VFE Layer-1 to describe the details in the following paragraph. Figure~\ref{fig:VFE} shows the architecture for VFE Layer-1.

\begin{figure}[!t]
\centering
    \includegraphics[width=0.8\linewidth]{./VFE_layer.pdf}
      \vspace{-0.3cm}
\caption{Voxel feature encoding layer.}
\label{fig:VFE}
\end{figure}

\begin{figure*}[!t]
\centering
    \includegraphics[width=0.8\linewidth]{RPN_new.pdf}
    \vspace{-0.3cm}
\caption{Region proposal network architecture.}
\label{fig:RPN}
\end{figure*}

Denote $\mathbf{V} = \{\mathbf{p}_i=[x_i, y_i, z_i, r_i]^{T} \in \mathbb{R}^4\}_{i=1 \ldots t}$ as a non-empty voxel containing $t \leq T$ LiDAR points, where $\mathbf{p}_i$ contains XYZ coordinates for the $i$-th point and $r_i$ is the received reflectance. We first compute the local mean as the centroid of all the points in $\mathbf{V}$, denoted as $(v_x, v_y, v_z)$. Then we augment each point $\mathbf{p}_i$ with the relative offset w.r.t. the centroid  and obtain the input feature set $\mathbf{V}_{\textrm{in}} = \{ \hat{\mathbf{p}}_i=[x_i, y_i, z_i, r_i, x_i-v_x, y_i-v_y, z_i-v_z]^{T} \in \mathbb{R}^7 \}_{i=1 \ldots t}$. Next, each $\hat{\mathbf{p}}_i$ is transformed through the fully connected network (FCN) into a feature space, where we can aggregate information from the point features $\mathbf{f}_i \in \mathbb{R}^m$ to encode the shape of the surface contained within the voxel. The FCN is composed of a linear layer, a batch normalization (BN) layer, and a rectified linear unit (ReLU) layer. After obtaining point-wise feature representations, we use element-wise MaxPooling across all $\mathbf{f}_i$ associated to $\mathbf{V}$ to get the locally aggregated feature $\tilde{\mathbf{f}} \in \mathbb{R}^m$ for $\mathbf{V}$. Finally, we augment each $\mathbf{f}_i$ with $\tilde{\mathbf{f}}$ to form the point-wise concatenated feature as $\mathbf{f}_i^{out} = [\mathbf{f}_i^{T}, \tilde{\mathbf{f}}^{T}]^{T} \in \mathbb{R}^{2m}$. Thus we obtain the output feature set $\mathbf{V}_{\textrm{out}} = \{ \mathbf{f}_i^{out} \}_{i \ldots t}$. All non-empty voxels are encoded in the same way and they share the same set of parameters in FCN.

We use $\textrm{VFE-i}(c_{in}, c_{out})$ to represent the $i$-th VFE layer that transforms input features of dimension $c_{in}$ into output features of dimension $c_{out}$. The linear layer learns a matrix of size $c_{in} \times (c_{out}/2)$, and the point-wise concatenation yields the output of dimension $c_{out}$.

Because the output feature combines both point-wise features and locally aggregated feature, stacking VFE layers encodes point interactions within a voxel and enables the final feature representation to learn descriptive shape information. The voxel-wise feature is obtained by transforming the output of VFE-$n$ into $\mathbb{R}^{C}$ via FCN and applying element-wise Maxpool where $C$ is the dimension of the voxel-wise feature, as shown in Figure~\ref{fig:VoxelNet_flowchart}.

\noindent{\bf Sparse Tensor Representation } By processing only the non-empty voxels, we obtain a list of voxel features, each uniquely associated to the spatial coordinates of a particular non-empty voxel. The obtained list of voxel-wise features can be represented as a sparse 4D tensor, of size $C \times D' \times H' \times W'$ as shown in Figure~\ref{fig:VoxelNet_flowchart}. Although the point cloud contains $\sim$100k points, more than $90\%$ of voxels typically are empty. Representing non-empty voxel features as a sparse tensor greatly reduces the memory usage and computation cost during backpropagation, and it is a critical step in our efficient implementation.

\subsubsection{Convolutional Middle Layers}
We use $\textrm{Conv}M\textrm{D}(c_{in}, c_{out}, \mathbf{k}, \mathbf{s}, \mathbf{p})$ to represent an $M$-dimensional convolution operator where $c_{in}$  and $c_{out}$ are the number of input and output channels, $\mathbf{k}$, $\mathbf{s}$, and $\mathbf{p}$ are the $M$-dimensional vectors corresponding to kernel size, stride size and padding size respectively. When the size across the $M$-dimensions are the same, we use a scalar to represent the size e.g. $k$ for $\mathbf{k} = (k,k,k)$.

Each convolutional middle layer applies 3D convolution, BN layer, and ReLU layer sequentially. The convolutional middle layers aggregate voxel-wise features within a progressively expanding receptive field, adding more context to the shape description. The detailed sizes of the filters in the convolutional middle layers are explained in Section~\ref{sec:training_details}.

\subsubsection{Region Proposal Network}

 Recently, region proposal networks~\cite{REF:NIPS2015_5638} have become an important building block of top-performing object detection frameworks~\cite{REF:DeepSlidingShapes,REF:cvpr17chen,REF:focalloss}. In this work, we make several key modifications to the  RPN architecture proposed in~\cite{REF:NIPS2015_5638}, and combine it with the feature learning network and convolutional middle layers to form an end-to-end trainable pipeline.

The input to our RPN is the feature map provided by the convolutional middle layers.  The architecture of this network is illustrated in Figure~\ref{fig:RPN}. The network has three blocks of fully convolutional layers. The first layer of each block downsamples the feature map by half via a convolution with a stride size of 2, followed by a sequence of convolutions of stride 1 ($\times q$ means $q$ applications of the filter). After each convolution layer, BN and ReLU operations are applied. We then upsample the output of every block to a fixed size and concatanate to construct the high resolution feature map. Finally, this feature map is mapped to the desired learning targets: (1) a probability score map and (2) a regression map.

\subsection{Loss Function}

Let $\{ a^{\text{\scriptsize pos}}_i \}_{i=1 \ldots N_{\text{\scriptsize pos}}}$ be the set of $N_{\text{\scriptsize pos}}$ positive anchors and $\{ a^{\text{\scriptsize neg}}_j \}_{j=1 \ldots N_{\text{\scriptsize neg}}}$ be the set of $N_{\text{\scriptsize neg}}$ negative anchors. We parameterize a 3D ground truth box as $(x_c^g, y_c^g, z_c^g, l^g, w^g, h^g, \theta^g)$, where $x_c^g, y_c^g, z_c^g$ represent the center location, $l^g, w^g, h^g$ are length, width, height of the box, and $\theta^g$ is the yaw rotation around Z-axis. To retrieve the ground truth box from a matching positive anchor parameterized as $(x_c^a, y_c^a, z_c^a, l^a, w^a, h^a, \theta^a)$, we define the residual vector $\mathbf{u}^* \in \mathbb{R}^7$  containing the 7 regression targets corresponding to center location $\Delta x, \Delta y, \Delta z$, three dimensions $\Delta l, \Delta w, \Delta h$, and the rotation $\Delta \theta$, which are computed as:
\begin{align}
\scriptsize
    &\Delta x = \frac{x_c^g - x_c^a}{d^a} , \Delta y = \frac{y_c^g - y_c^a}{d^a} ,
    \Delta z = \frac{z_c^g - z_c^a}{h^a} , \nonumber \\
    &\Delta l = \log(\frac{l^g}{l^a}) , \Delta w = \log(\frac{w^g}{w^a}) , \Delta h = \log(\frac{h^g}{h^a}) , \\
    &\Delta \theta = \theta^g - \theta^a \nonumber
\label{eqn:reg_targets}
\end{align}
where $d^a = \sqrt{(l^a)^2 + (w^a)^2}$ is the diagonal of the base of the anchor box. Here, we aim to directly estimate the oriented 3D box and normalize $\Delta x$ and $\Delta y$ homogeneously with the diagonal $d^a$, which is different from \cite{REF:NIPS2015_5638,REF:DeepSlidingShapes,REF:VeloFCN, REF:3DFCN,REF:nips15chen,REF:cvpr16chen,REF:cvpr17chen}. We define the loss function as follows:

\begin{eqnarray}
    L & = & \alpha\frac{1}{N_{\text{\scriptsize pos}}} \sum\limits_{i} L_{\text{\scriptsize cls}}(p^{\text{\scriptsize pos}}_i, 1) +
            \beta\frac{1}{N_{\text{\scriptsize neg}}} \sum\limits_{j} L_{\text{\scriptsize cls}}(p^{\text{\scriptsize neg}}_j, 0) \nonumber \\
            & + & \frac{1}{N_{\text{\scriptsize pos}}} \sum\limits_{i} L_{\text{\scriptsize reg}}(\mathbf{u}_i, \mathbf{u}_i^*)
\label{eqn:loss_function}
\end{eqnarray}
where $p^{\text{\scriptsize pos}}_i$ and $p^{\text{\scriptsize neg}}_j$ represent the softmax output for positive anchor $a^{\text{\scriptsize pos}}_i$ and negative anchor $a^{\text{\scriptsize neg}}_j$ respectively, while $\mathbf{u}_i \in \mathbb{R}^7$ and $\mathbf{u}_i^* \in \mathbb{R}^7$ are the regression output and ground truth for positive anchor $a^{\text{\scriptsize pos}}_i$. The first two terms are the normalized classification loss for $\{ a^{\text{\scriptsize pos}}_i \}_{i=1 \ldots N_{\text{\scriptsize pos}}}$ and $\{ a^{\text{\scriptsize neg}}_j \}_{j=1 \ldots N_{\text{\scriptsize neg}}}$, where the $L_{\text{\scriptsize cls}}$ stands for binary cross entropy loss and $\alpha$, $\beta$ are postive constants balancing the relative importance. The last term $L_{reg}$ is the regression loss, where we use the SmoothL1 function~\cite{REF:Girshick:2015,REF:NIPS2015_5638}.

\begin{figure}[!t]
\centering
    \includegraphics[width=0.75\linewidth]{Implementation.pdf}
\caption{Illustration of efficient implementation.}
\vspace{-0.3cm}
\label{fig:Implementation}
\end{figure}

\subsection{Efficient Implementation}
\label{subsec:implementation}

GPUs are optimized for processing dense tensor structures. The problem with working directly with the point cloud is that the points are sparsely distributed across space and each voxel has a variable number of points.  We devised a method that converts the point cloud into a dense tensor structure where stacked VFE operations can be processed in parallel  across points and voxels.

The method is summarized in Figure~\ref{fig:Implementation}. We initialize a $K \times T \times 7$ dimensional tensor structure to store the voxel input feature buffer where $K$ is the maximum number of non-empty voxels, $T$ is the maximum number of points per voxel, and $7$ is the input encoding dimension for each point. The points are randomized before processing. For each point in the point cloud, we check if the corresponding voxel already exists. This lookup operation is done efficiently in $O(1)$ using a hash table where the voxel coordinate is used as the hash key. If the voxel is already initialized we insert the point to the voxel location if there are less than $T$ points, otherwise the point is ignored. If the voxel is not initialized, we initialize a new voxel, store its coordinate in the voxel coordinate buffer, and insert the point to this voxel location. The voxel input feature and coordinate buffers can be constructed via a single pass over the point list, therefore its complexity is $O(n)$.  To further improve the memory/compute efficiency it is possible to only store a limited number of voxels ($K$) and ignore points coming from voxels with few points.

After the voxel input buffer is constructed, the stacked VFE only involves point level and voxel level dense operations which can be computed on a GPU in parallel. Note that,  after concatenation operations in VFE, we reset the features corresponding to empty points to zero such that they do not affect  the computed voxel features. Finally, using the stored coordinate buffer we reorganize the computed sparse voxel-wise structures to the dense voxel grid. The following  convolutional middle layers and RPN operations work on a dense voxel grid which can be efficiently implemented on a GPU.


\section{Introduction}
NASA's {\em Kepler Space Telescope}, operating in its prime mission (2009--2013; \citealt{Borucki10a}), shed light on many fundamental properties of exoplanets. Among these are the occurrence of planets as small as Earth around Sun-like and low-mass stars (e.g. \citealt{Petigura13b,Dressing15}) and the diversity of planetary bulk compositions \citep{Marcy14,Weiss14,Rogers15} extending down to Earth-size (e.g. \citealt{Howard13,Jontof-Hutter15}).

Now operating in its two-wheel \textit{K2}\xspace mode \citep{Howell14}, \Kepler observes a different region of sky every three months. \textit{K2}\xspace is conducting a wider, more shallow survey that complements the narrow, deep survey of the prime mission. Among its many accomplishments to date, \textit{K2}\xspace has significantly increased the number of transiting planets around moderately bright stars \citep{Crossfield16}, which will enable more detailed studies of exoplanet bulk composition using precision radial velocity facilities. \textit{K2}\xspace has also revealed planets around newborn stars \citep{David16,Mann16} and planets around white-dwarfs \citep{Vanderburg15WD}. Due to community-driven target selection, a large fraction of the \textit{K2}\xspace targets are M-dwarfs, resulting in the detection of planets in or near the habitable zone (e.g. \citealt{Crossfield15,Montet15,Petigura15,Schlieder16}).

In this paper, we provide a catalog of planet candidates and eclipsing binaries from the second year of \textit{K2}\xspace operations, corresponding to campaigns 5--8 (C5--C8). Section~\ref{sec:planet-candidate-search} presents our methodology for correcting spacecraft systematics in \textit{K2}\xspace photometry and identifying planet candidates. In Section~\ref{sec:spectroscopy}, we describe our spectroscopic follow-up program and present refined stellar parameters enabled by these spectra. We present our catalog of planet candidates and eclipsing binaries in Section~\ref{sec:planet-candidate-table} and summarize our findings in Section~\ref{sec:conclusions}.

\section{Identifying Planet Candidates}
\label{sec:planet-candidate-search}

\subsection{Photometry}
\label{sec:photometry}
During its prime mission, \Kepler achieved photometric precisions of $\approx$40~ppm on 6.5 hour timescales \citep{Christiansen12} for targets of $\approx$12~mag in the \Kepler bandpass
(i.e. $\textit{Kp}\xspace\approx12$~mag). For many stars, photometric precision was limited by intrinsic stellar variability rather than photon-limited or instrumental errors. This exquisite precision was due in large part to stable pointing enabled by four (and later three) reaction wheels which stabilized the telescope against solar radiation pressure across the three axes of the telescope. Photometry for individual target stars was extracted using stationary software apertures composed of integer numbers of connected pixels.

During \textit{K2}\xspace operations, where the spacecraft uses the two remaining operational reaction wheels, solar radiation pressure causes drifts of $\sim$1~pixel to occur on $\sim$6~hr timescales. As stars drift across the CCD, variations in pixel sensitivities and variable aperture losses result in apparent brightness variations.

Several techniques have been developed to correct for the position-dependent brightness variations due to the unstable platform of \textit{K2}\xspace. A non-exhaustive list includes {\tt k2sff} \citep{Vanderburg14}, {\tt k2phot} \citep{Crossfield15,Petigura15,Crossfield16}, and {\tt k2sc} \citep{Aigrain15} which model stellar brightness as a function of spacecraft orientation with a function and remove it from the light curve. The {\tt everest} package \citep{Luger16} builds on the pixel-level decorrelation (PLD) approach developed for \Spitzer \citep{Deming15} and decorrelates against the pixel-by-pixel photometric timeseries.

We generated light curves for \stat{n-lightcurves} stars observed by \textit{K2}\xspace during C5--C8 using the publicly-available {\tt k2phot} Python package.%
\footnote{https://github.com/petigura/k2phot (commit a0d507)}
The general methodology is described in previous works \citep{Crossfield15,Petigura15,Crossfield16}. However, due to the evolving nature of \textit{K2}\xspace systematics, we have continued to adapt and refine {\tt k2phot} and summarize the methodology below.

Since systematics in the photometry are largely due to pointing drifts, accurate knowledge of the spacecraft orientation is important for removing orientation-dependent systematics. We characterize the time-dependent orientation of the \Kepler spacecraft by analyzing the positions of $\approx$100 bright but unsaturated stars having $\ensuremath{Kp}\xspace\approx12$~mag on a representative output channel of the \Kepler CCD.%
\footnote{The \Kepler CCD contains 84 output channels (Kepler Instrument Handbook; Van Cleve et al. 2016), of which 76 were operational during C5-C8. An additional module (4 output channels) failed during C10.}
For each long-cadence integration, we solve for the affine transformation that maps that frame to an arbitrary reference frame. We then use the sequence of affine transformations to transform a reference pixel coordinate on a reference frame%
\footnote{For example, a pixel on row 500, column 500, and frame 2000.}
to the pixel coordinate on all other frames in the campaign.

We extract photometry using stationary apertures. Aperture size is described by a single variable \ensuremath{N_{\mathrm{pix}}}\xspace, the number of pixels in the aperture. Apertures are constructed to accommodate image motion during a campaign. We construct apertures using a composite image constructed from the 90th percentile intensity value of all frames in a campaign. Because the stars move during \textit{K2}\xspace observations, the 90th percentile image is smeared out and the apertures constructed from this image accommodate the drifts, mitigating severe aperture losses. The apertures are then constructed by selecting the pixel closest to the expected position of the target star, as predicted by the WCS coordinates provided by the \textit{K2}\xspace project. Additional pixels are added iteratively by selecting the brightest pixel touching the current mask.

During the photometric extraction, we search for the aperture size $N_\mathrm{pix,min}$ that minimizes noise on three-hour times scales. This aperture size strikes a balance between the desire to minimize systematic noise which grows with {\em decreasing} aperture size, and background noise, which grows with {\em increasing} aperture size. We select an initial size $N_\mathrm{pix,0}$, which is motivated by previous analyses of stars with similar \ensuremath{Kp}\xspace. We then try six logarithmically-spaced \ensuremath{N_{\mathrm{pix}}}\xspace between $N_\mathrm{pix,0} / 4$ and $N_\mathrm{pix,0} \times 4$, which samples a curve describing noise as a function of \ensuremath{N_{\mathrm{pix}}}\xspace. We find $N_\mathrm{pix,min}$ using upto three iterations of Newton's method. While testing different aperture sizes, we constrain  \ensuremath{N_{\mathrm{pix}}}\xspace to be between nine pixels and the total number of pixels in the target pixel file.

After extraction of the photometry we have a sequence of flux as a function of $x$, $y$, and $t$. We model out changes in flux that correlate with changes in $x$, $y$, and $t$ using a Gaussian process with a squared-exponential covariance kernel, which is characterized by the following seven hyper-parameters $A_x$, $l_x$, $A_y$, $l_y$, $A_t$, $l_t$, and $\sigma$. Here, $A$ corresponds to the amplitude of the GP, $l$ corresponds to characteristic length scale, and $\sigma$ accounts for a white noise component. Choosing the appropriate hyper-parameters can be computationally intensive on a star-by-star basis. We therefore adopt a scheme from \cite{Aigrain15} which optimizes the hyper-parameters on a subset of the photometry using a differential evolution global optimizer \citep{Storn97}.

We produced light curves for \stat{n-lightcurves} stars in C5--C8, which are available on the Exoplanet Follow-up Observing Program (ExoFOP) website.%
\footnote{https://exofop.ipac.caltech.edu/}
Along with the photometry, we included photometric diagnostic plots showing the extraction aperture and resulting detrended light curve. Figure~\ref{fig:211736671} shows these diagnostic plots for an example planet candidate around EPIC-211736671.

\begin{figure*}
\gridline{\fig{fig_211736671-aper.pdf}{0.9\textwidth}{(a)}}
\gridline{\fig{fig_211736671-photometry.pdf}{0.9\textwidth}{(b)}}
\caption{An example of the photometric diagnostic plots included as standard data products on the ExoFOP. Panel {\bf (a)}: three images of an example planet candidate (EPIC-211736671), observed in C5. {\em Left:} the median of all long cadence C5 frames with the optimal extraction aperture shown in green. {\em Center:} the same region of sky as observed by the first Palomar Observatory Sky Survey (POSS-I). The orange region corresponds to the boundaries of the K2 frame. {\em Right:} the same region of sky as observed by POSS-II. Panel {\bf (b):} aperture photometry before and after subtraction of our systematic noise model.\label{fig:211736671}}

\end{figure*}

\subsection{Transiting Planet Search}
Our general transit search and vetting process is similar to that described in \cite{Crossfield16}. We give a brief summary of here. We searched the calibrated photometry for transiting planets using the publicly-available \texttt{TERRA}\xspace algorithm.%
\footnote{https://github.com/petigura/terra (commit 9739e9)}
\texttt{TERRA}\xspace is a matched filter-based approach and is described by \cite{Petigura13b}. \texttt{TERRA}\xspace convolves the photometry with a box-shaped approximation of a transit profile to compute a Single Event Statistic (SES) at every \textit{K2}\xspace long-cadence measurement. The SES time series is phase-folded according to a finely-spaced grid of trial periods $P$ and times of first transit \T{0}.

A classical matched-filter algorithm would then compute a Multiple Event Statistic (MES) by summing SES at each trial (\per, \T{0}), which is optimal given uncorrelated Gaussian noise . However, in \textit{K2}\xspace photometry we observe more frequent non-Gaussian anomalies relative to \Kepler prime photometry due to the aggressive detrending that must be performed in order to remove the instrumental systematics described above. A traditional MES computation resulted in many spurious peaks with apparently high MES, but were later easily identifiable as anomalies through inspection. We address outliers by calculating MES only after removing the two highest SES peaks at each trial (\per, \T{0}).  Spurious peaks due to the chance alignment of two outliers are eliminated. This non-linear filter removes many spurious detections and eases the burden during manual vetting, described below. One consequence is that \texttt{TERRA}\xspace does not identify planets with one or two transits occurring in a \textit{K2}\xspace campaign. Such transits are sometimes identified by eye, but we caution that many are likely overlooked. These events are especially amenable to searches by citizen scientists (see e.g., Christiansen et al. {\em accepted in AJ}).%
\footnote{\url{https://www.zooniverse.org/projects/ianc2/exoplanet-explorers}}

\texttt{TERRA}\xspace identifies $\sim$1000 Threshold-Crossing Events (TCEs) per campaign. A TCE is a particular combination of (\per, \T{0}) that has MES that exceeds some threshold, but may not be an astrophysical transit. If a candidate has a periodic dimming of a consistent shape it is elevated to the status of a ``\textit{K2}\xspace Object of Interest'' (K2OIs) which are consistent with an astrophysical transit or eclipse.

Our team visually inspects each K2OI to look for a robust indication that the target is an eclipsing binary. We look for secondary eclipses, which indicate that the transiting object is self-luminous and not a planet. Secondary eclipses associated with binaries with circular orbits are shifted in phase from the primary eclipse by $180^{\circ}$. We search for secondary eclipses at all phases to allow for eccentric orbits. We also look for obvious odd/even variations, which indicate that \texttt{TERRA}\xspace has identified a nearly circular EB at half the orbital period. We also identify stars that show variability that is phase-locked to the eclipse, which is a strong indicator of star-star modulation (ellipsoidal, reflection, or relativistic beaming).

Our eclipsing binary designation does not incorporate transit depth or whether the light curve is V-shaped. While these attributes are strong indicators of EB status, they are not conclusive. Planets transiting small M-dwarf stars can easily produce transits deeper than 1\% and short period transits may appear V-shaped due to the 30-minute sampling of \textit{K2}\xspace. We defer a detailed false positive calculation for a later paper.

In total, we identified \stat{n-k2oi} K2OIs, associated with $\stat{n-k2oi-stars}$ stars. Of these, \stat{n-k2oi-eb}/\stat{n-k2oi} are likely eclipsing binaries, and we refer to the remaining \stat{n-k2oi-not-eb} as planet candidates.

\subsection{Light Curve Fitting}
We fit the calibrated photometry according to the methodology of \cite{Crossfield16}. In brief, we used the publicly-available {\tt batman} light curve code \citep{Kreidberg15} to generate model light curves which we then compared against the photometry. We first derived a maximum likelihood solution and then derived parameter uncertainties using Markov Chain Monte Carlo (MCMC).%
\footnote{using the affine-invariant sampler of \cite{Goodman10} as implemented in Python by \cite{Foreman-Mackey15}}

In our modeling, the following parameters are allowed to vary: time of first transit $T_0$, orbital period $P$, inclination $i$, scaled semi-major axis $a/\Rstar$, planet-star radius ratio \ensuremath{R_P}\xspace/\Rstar, orbital eccentricity $e$, longitude of periastron $\omega$, linear limb-darkening coefficient $u$, fractional light curve dilution $\delta$, and the out-of-transit flux level.

During the fitting, we adopted the following priors:

\begin{itemize}
\item {\em Period}. Gaussian prior centered on maximum likelihood $P$ having dispersion of 0.01~days.%
\footnote{We imposed weak priors on $P$ and $T_0$ to keep the MCMC walkers from jumping too far from the likelihood mode and wandering away. After performing the MCMC exploration, we verified that $P$ and $T_0$ were more tightly constrained by the photometry than by the priors. The uncertainties on $T_0$ are typically 2.5\% of the prior width (median value) and no more than 50\% the prior width. The uncertainties on P are typically 2\% of the prior width (median value) and no more than 60\% the prior width.}
\item {\em Time of transit}. Uniform prior centered on maximum likelihood $T_0$ having dispersion of $0.06 \times P$.
\item {\em Radius ratio}. Uniform prior, $\ensuremath{R_P}\xspace/\Rstar$ = [$-1$,$+1$]. Following \cite{Eastman13}, we allow for negative $\ensuremath{R_P}\xspace/\Rstar$ in our sampling to avoid the Lucy-Sweeney-type bias that results from treating $\ensuremath{R_P}\xspace/\Rstar$ as a positive-definite quantity \citep{Lucy71}.
\item {\em Eccentricity}. Gaussian prior centered at $10^{-4}$ having dispersion of $10^{-3}$. This effectively restricts the orbits to circular.
\item {\em Longitude of periastron}. Uniform prior, $\omega$ = [0,$2\pi$].
\item {\em Inclination}. Uniform prior, $i$ = [$50^{\circ}$,$90^{\circ}$].
\item {\em Limb-darkening}. Gaussian prior on $u$ where the mean and dispersion are computed using the publicly-available Limb-Darkening Toolkit (LDTK; \citealt{Parviainen15}). LDTK computes the distribution of $u$ given Gaussian constraints on \teff, \logg, and \ensuremath{\mathrm{[Fe/H]}}\xspace. For consistency, we used \teff, \logg, and \ensuremath{\mathrm{[Fe/H]}}\xspace constrained by broadband photometry from \cite{Huber16}. While photometrically-constrained \logg and \ensuremath{\mathrm{[Fe/H]}}\xspace are low-precision, $u$ is only weakly dependent on these parameters, and the derived transit parameters are only weakly dependent on $u$.
\item {\em Dilution}. Log-uniform prior, $\log \delta$ = [$10^{-6}$,$10^{0}$]. Our fits do not incorporate dilution constraints, so $\delta$ always reverts to the prior. We include $\delta$ so we can incorporate dilution constraints at later times.%
\footnote{The non-zero prior on $\delta$ slightly alters the derived value of $\ensuremath{R_P}\xspace/\Rstar$. However, because the median $\delta$ is $10^{-3}$, this amounts to a fractional change in derived radius ratio of $\Delta(\ensuremath{R_P}\xspace/\Rstar)/ (\ensuremath{R_P}\xspace/\Rstar) \approx \delta/ 2 = 5\times10^{-4}$, which may be ignored.}
\end{itemize}

In Table~\ref{tab:candidates}, we report 1$\sigma$ credible ranges on $P$, $T_0$, $\ensuremath{R_P}\xspace/\Rstar$, transit duration $T_{14}$, and impact parameter $b$.

\section{Spectroscopy}
\label{sec:spectroscopy}
\subsection{Spectroscopic Follow Program}
As a part of our team's standard follow-up efforts, we obtained optical spectra of $\stat{n-hires-spec}$ C5--C8 target stars using the HIgh Resolution Echelle Spectrometer (HIRES; \citealt{Vogt94}) on the Keck-I 10~m telescope. We gathered spectra for the purpose of improving host star parameters and to place limits on the presence of companion stars with small separations through searches for spectroscopic binaries (SB2s). We aimed to obtain a spectrum of every K2OI brighter than $V = 14.0$~mag. For G stars, this limit corresponds roughly to $\ensuremath{Kp}\xspace = 13.6$~mag.

Table~\ref{tab:spec} lists the C5--C8 targets that we observed with HIRES, along with the results from our stellar characterization and search for spectroscopic binaries, which are described in Sections~\ref{sec:stellar-characterization} and \ref{sec:sb2}, respectively. We obtained HIRES spectra of \stat{n-hires-spec-pc}/\stat{n-k2oi-not-eb-stars} of the planet candidate host stars and for \stat{n-hires-spec-eb}/\stat{n-k2oi-eb} of the likely EBs. In addition, we observed \stat{n-hires-spec-other} other C5--C8 targets that we did not identify as candidates. These were typically observed because they were identified as planet candidates by other groups.

We used the HIRES exposure meter to obtain consistent SNR levels depending on $V$-band apparent magnitude: SNR = 45 ($V < 13.0$~mag), SNR = 32 ($V = 13.0-14.0$~mag), and SNR = 20 ($V > 14.0$~mag). Exposure times ranged from $\approx$10~s for $V$ = 9~mag targets to $\approx$400~s for $V$ = 15~mag targets. SNR is computed per reduced HIRES pixel on blaze at 5500~\AA. Our HIRES follow-up was nearly complete down to $V = 14$~mag. Figure~\ref{fig:kepmag-color} shows the distribution of candidate hosts as a function of \ensuremath{Kp}\xspace and \jmk color. The candidates with HIRES spectra are labeled. Figure~\ref{fig:smemp-lincomb} shows a spectral segment for one K2OI to illustrate typical spectral resolution and SNR level.

\begin{figure}
\centering
\includegraphics[width=0.5\textwidth]{fig_kepmag-color.pdf}
\caption{Distribution of K2OIs as a function of \ensuremath{Kp}\xspace and \jmk color. The colored/gray points represent targets with/without a HIRES spectrum. The marker color and shape represent a target's specific \textit{K2}\xspace campaign. The HIRES follow-up is nearly complete to \textit{Kp}\xspace~=~14~mag.\label{fig:kepmag-color}}
\end{figure}

\begin{figure*}
\centering
\includegraphics[width=0.9\textwidth]{fig_EPIC212006344-lincomb.pdf}
\caption{Example HIRES characterization spectrum with \texttt{SpecMatch-Emp}\xspace fit. Blue spectrum is of EPIC212006344, an M0 dwarf, and illustrates the typical SNR from our characterization program of $\approx45$/pixel. The spectrum contains a dense forest of molecular lines, making {\em ab initio} spectral synthesis challenging. Red spectra are drawn from \texttt{SpecMatch-Emp}\xspace library and are identified as similar to the target spectrum by the \texttt{SpecMatch-Emp}\xspace algorithm. The green and black spectra are the best-fitting linear combination spectrum and residual spectrum, respectively.\label{fig:smemp-lincomb}}
\end{figure*}

\subsection{Stellar Characterization}
\label{sec:stellar-characterization}
We used our spectra to improve the precision of stellar and planetary properties such as \Rstar and \ensuremath{R_P}\xspace. We analyzed each spectrum with one of two related publicly-available codes: \texttt{SpecMatch-Syn}\xspace \citep{Petigura15thesis}%
\footnote{\url{https://github.com/petigura/specmatch-syn}}
and \texttt{SpecMatch-Emp}\xspace \citep{Yee17}.%
\footnote{\url{https://github.com/samuelyeewl/specmatch-emp}}

\texttt{SpecMatch-Syn}\xspace fits five regions of optical spectrum by interpolating within a grid of model spectra from \cite{Coelho05}. Recently, \texttt{SpecMatch-Syn}\xspace enabled a homogeneous analysis of 1305 spectra of planet hosts identified during the prime \Kepler mission \citep{Petigura17b}. For stars with \teff = 4700--6500~K and \vsini < 20~km s$^{-1}$\xspace, \texttt{SpecMatch-Syn}\xspace achieves precisions of 60~K in \teff, 0.10~dex in \logg, and 0.04~dex in \ensuremath{\mathrm{[Fe/H]}}\xspace, and 1~km s$^{-1}$\xspace in \vsini.

We converted \teff, \logg, and \ensuremath{\mathrm{[Fe/H]}}\xspace into \ensuremath{M_{\star}}\xspace and \Rstar using the publicly-available \texttt{isoclassify}\xspace Python package \citep{Huber17},%
\footnote{\url{https://github.com/danxhuber/isoclassify}}
which uses the MESA Isochrones and Stellar Tracks (MIST) database \citep{Choi16,Paxton11,Paxton13,Paxton15}. While \texttt{SpecMatch-Syn}\xspace returns \teff precise to 60~K, as tested against other spectral synthesis codes, there are known offsets between spectroscopic \teff and other techniques such as the Infrared Flux Method (IRFM) and interferometry. For a detailed discussion of different \teff scales, see \cite{Brewer16}. To account for systematic uncertainties associated with the spectroscopic \teff scale, we have increased the \teff uncertainties to 100~K during the isochrone modeling.

The radius uncertainties derived from \texttt{SpecMatch-Syn}\xspace parameters do not incorporate uncertainties associated with the MIST models themselves. \cite{Johnson17} estimated the size of these model-dependent uncertainties through a comparison of stellar radii derived using MIST models and Dartmouth Stellar Evolution Program (DSEP) models \citep{Dotter08} with identical inputs. They estimated that model-dependent radius errors are $\approx$2\% for dwarf stars ($\logg < 3.9$) and $\approx$10\% for evolved stars ($\logg>3.9$). These model-dependent uncertainties are typically smaller than the formal radius uncertainties returned by \texttt{isoclassify}\xspace.

For $\teff \lesssim 4700$~K, \texttt{SpecMatch-Syn}\xspace does not return reliable parameters, due to the onset of molecular lines that are not well-treated in the \cite{Coelho05} models. While the high-resolution optical spectra of stars later than $\sim$K4 is challenging to compute directly, their spectra contain a wealth of information, which can be used to constrain stellar properties. \texttt{SpecMatch-Emp}\xspace circumvents the challenges in spectral synthesis by matching a target spectrum against an empirical spectra library of $\approx$400 touchstone stars with well-known parameters measured through other methods such as SED-fitting, interferometry, or the IRFM. \texttt{SpecMatch-Emp}\xspace interpolates between this library of empirical spectra to find linear combination of library spectra that best reproduces the target spectrum. \texttt{SpecMatch-Emp}\xspace achieves precisions of 70~K in \teff, 10\% in \Rstar, and 0.12~dex in \ensuremath{\mathrm{[Fe/H]}}\xspace.

We adopt parameters from \texttt{SpecMatch-Syn}\xspace for stars hotter than 4700~K%
\footnote{as measured by \texttt{SpecMatch-Emp}\xspace}
and \texttt{SpecMatch-Emp}\xspace for cooler stars. Figure~\ref{fig:spec-hr-diagram} shows the \teff and \Rstar for K2OIs with reliable spectroscopic parameters. Our adopted stellar parameters are listed in Table~\ref{tab:spec}. Our team also conducts a parallel characterization of cool stars using NIR spectroscopy. Stellar properties up through campaign 7 are given in \cite{Martinez17} and \cite{Dressing17}.

\subsection{Searches for Spectroscopic Binaries}
\label{sec:sb2}

Each HIRES spectrum is methodically searched for secondary spectral lines using the \texttt{ReaMatch}\xspace algorithm \citep{Kolbl15}.  To identify secondary spectra, each spectrum is first cross-correlated against a set of previously-observed spectra. This catalog has spectra with \teff = 3500--6500~K and \logg = 3.0--4.5~dex. The best-matching spectrum is subtracted from the target spectrum and the residuals are again cross-correlated with the catalog spectra. \texttt{ReaMatch}\xspace is sensitive to companions down to 1\% the brightness of the primary having RV offsets $|\Delta v|$~>~10~km s$^{-1}$\xspace. Optimized for slowly-rotating FGKM stars, \texttt{ReaMatch}\xspace is insensitive to SB2s orbiting primaries with \vsini > 10 km s$^{-1}$\xspace or with \teff outside 3500--6500 K. Table~\ref{tab:spec} lists the results of our SB2 search.

\begin{figure*}
\centering
\includegraphics[width=0.75\textwidth]{fig_EPIC212006344-library.pdf}
\caption{\texttt{SpecMatch-Emp}\xspace characterization of EPIC212006344. Black points: \teff, \Rstar, and \ensuremath{\mathrm{[Fe/H]}}\xspace from the \texttt{SpecMatch-Emp}\xspace library. Green triangles: properties of the closest-matching library spectra (red spectra in Figure~\ref{fig:smemp-lincomb}). Purple square: linear combination of library spectra spectra. The final derived parameters are \teff~=~$3925\pm70$~K, \ensuremath{\mathrm{[Fe/H]}}\xspace = $0.43\pm0.12$~dex, \Rstar~=~$0.63\pm0.10$~dex.\label{fig:smemp-library}}
\end{figure*}

\begin{figure}
\centering
\includegraphics[width=0.48\textwidth]{fig_spec-hr-diagram.pdf}
\caption{K2OIs with reliable effective temperatures and stellar radii, as measured by \texttt{SpecMatch-Syn}\xspace (blue) or \texttt{SpecMatch-Emp}\xspace (green). Median uncertainties are shown at top right.\label{fig:spec-hr-diagram}}
\end{figure}

\section{Planet Candidates}
\label{sec:planet-candidate-table}

We list the \stat{n-k2oi-not-eb} planet candidates and \stat{n-k2oi-eb} likely EBs in Table~\ref{tab:candidates}. We compute planetary radius by combining \rrat measured from the light curve with the best available \Rstar. For the stars without spectra, we estimate \Rstar from broadband photometry, according to the following procedure: we estimate spectral types (SpTs) using tabulated photometric relations \citep{Kraus07,Pecaut13,Rodriguez13} and convert SpTs into \Rstar based on interferometric studies \citep{Boyajian12b}. These stellar radii are crude and we estimate their typical uncertainties to be $\approx$40\%, typical of errors derived from broadband photometry \citep{Brown11}.

Figure~\ref{fig:per-pradius} shows the distribution of the \stat{n-k2oi-not-eb} planet candidates in the $P-\ensuremath{R_P}\xspace$ plane. Figure~\ref{fig:hist-per-prad-kepmag} shows 1D histograms of our candidates as a function of $P$, $\ensuremath{R_P}\xspace$, and \ensuremath{Kp}\xspace. The median host star is nearly two magnitudes brighter in the \Kepler bandpass than the median KOI from the prime \Kepler mission (\stat{kepmag-k2oi-not-eb}~mag vs \stat{kepmag-koi}~mag, \citealt{Mullally15}). Our candidates have the following multiplicity distribution: \stat{n-k2oi-not-eb-stars-1} singles, \stat{n-k2oi-not-eb-stars-2} doubles, and \stat{n-k2oi-not-eb-stars-3} triple planet system.

We consulted the NASA Exoplanet Archive (NEA; \citealt{Akeson13})%
\footnote{https://exoplanetarchive.ipac.caltech.edu/}
to check whether previous analyses have reported significant numbers of candidates presented in this work. Of the catalogs incorporated into the NEA as of 2017-11-09, \citeauthor{Barros16} (2016; B16 hereafter) and  \citeauthor{Pope16} (2016; P16 hereafter) included 10 or more candidates from C5--C8.

B16 reported 172 planet candidates from C1--C6, of which 86 were in C5 and C6. Our catalog contains 107 candidates from C5 and C6. The two catalogs share 49 candidates. There are 58 candidates in our catalog that are not in B16, and there are 37 candidates in B16 that are not in our catalog.

P16 reported 168 planet candidates in C5 and C6. Of these, our catalog includes 73 candidates. There are 34 candidates in our catalog that are not in P16, and there are 95 candidates in P16 that are not in our catalog.

As a final point of comparison, B16 and P16 share 59 candidates. There are 27 candidates in B16 that are not in P16, and there are 109 candidates in P16 that were not in B16. Figure~\ref{fig:venn} is a Venn diagram that summarizes the degree of overlap between the various samples.

A detailed analysis of why any particular candidate appeared in one catalog and not another is beyond the scope of this work. Broadly-speaking, the lack perfect overlap likely arises due to differences in photometric extraction algorithms, transit search algorithms, adopted signal-to-noise threshold for candidate status, and vetting diagnostics.

\begin{figure}
\centering
\includegraphics[width=0.5\textwidth]{fig_per-pradius.pdf}
\caption{Planet size and orbital period for $\stat{n-k2oi-not-eb}$ planet candidates identified by \textit{K2}\xspace in C5--C8. We have excluded targets identified as likely EBs from their light curve morphology. The legend at lower right links marker shape/color to a specific campaign and gives the total number of candidates identified in each campaign.\label{fig:per-pradius}}
\end{figure}

\begin{figure*}[ht!]
\begin{center}
\includegraphics[trim=2.cm 0.5cm 0.5cm 0cm,width=0.32\textwidth]{fig_hist-per.pdf}
\hspace{0.1cm}
\includegraphics[trim=0.5cm 0.5cm 1.5cm 0cm,width=0.32\textwidth]{fig_hist-prad.pdf}
\includegraphics[trim=0.5cm 0.5cm 1.5cm 0cm,width=0.32\textwidth]{fig_hist-kepmag.pdf}
\caption{{\em Left:} Orbital periods of transit-like signals identified in our analysis. The orange histogram (axis at left) indicates the distribution of planet candidates. The pale, narrow-binned histogram (axis at right) indicates the Threshold-Crossing Events (TCEs) identified by \texttt{TERRA}\xspace in our initial transit search with MES\,$\ge10$.
{\em Middle:} Cumulative histograms of radii for our planet candidates. Most are moderately bright at $Kp \le 13$~mag, but at large radii over half orbit fainter stars; a large fraction of this second group are likely false positives.
{\em Right:} The orange histogram (axis at left) shows the distribution of \ensuremath{Kp}\xspace for planet candidates. For comparison, the pale histogram (axis at right) shows all target stars from C5--C8.%
\label{fig:hist-per-prad-kepmag}}
\end{center}
\end{figure*}

\begin{figure*}[ht!]
\centering
\includegraphics[trim=1.75cm 1cm 1.75cm 0cm, width=0.32\textwidth]{fig_venn-p17-b16.pdf}
\includegraphics[trim=1.75cm 1cm 1.75cm 0cm,width=0.32\textwidth]{fig_venn-p17-p16.pdf}
\includegraphics[trim=1.75cm 1cm 1.75cm 0cm,width=0.32\textwidth]{fig_venn-b16-p16.pdf}
\caption{Venn diagrams that compare the candidates from C5 and C6 reported in the following catalogs: This work, \cite{Barros16}, and \cite{Pope16}.\label{fig:venn}}
\end{figure*}

\section{Summary}
\label{sec:conclusions}

We report \stat{n-k2oi-not-eb} planet candidates orbiting \stat{n-k2oi-not-eb-stars} stars detected through a systematic search of \textit{K2}\xspace photometry from campaigns C5--C8. We also identified \stat{n-k2oi-eb} likely EBs based on their light curve morphology. We obtained Keck/HIRES optical spectra of \stat{n-hires-spec-pc}/\stat{n-k2oi-not-eb-stars} planet candidate host stars and \stat{n-hires-spec-eb}/\stat{n-k2oi-eb} EBs to improve our understanding of host star and planet properties and to search for binary companions.

A small fraction of our planet candidates reside in multi-candidate systems (\stat{n-k2oi-not-eb-stars-2} doubles and \stat{n-k2oi-not-eb-stars-3} triple). These systems have a low false positive probability ($\lesssim 1\%$) due to their multiplicity \citep{Lissauer12}. The remaining \stat{n-k2oi-not-eb-stars-1} candidates are well-vetted and well-characterized planet candidates, but have yet been confirmed or statistically validated. Statistical validation requires a detailed analysis of light curve shape and constraints on the presence of blends from high-resolution imaging. \cite{Crossfield16} performed such an analysis to validate 104 K2 planet candidates identified during C0--C4. Our team's high contrast imaging followup will be presented in Gonzales et al. {\em in prep.} An analysis of the false positive probabilities of our candidates will be presented in Livingston et al. {\em in prep.}

Our typical candidate is two magnitudes brighter than the typical candidate from the \Kepler prime mission due to the larger region of sky observed by \textit{K2}\xspace. As a result, these candidates make up a valuable sample for further characterization efforts, such as Doppler measurements of planet masses.

\software{batman \citep{Kreidberg15}, SpecMatch-Syn \citep{Petigura15thesis}, SpecMatch-Emp \citep{Yee17}, k2phot (https://github.com/petigura/k2phot), isoclassify \citep{Huber17}, isochrones \citep{Morton15iso}.}

\acknowledgments
We thank the anonymous referee for a thoughtful reading of the manuscript and for useful suggestions. EAP acknowledges support from Hubble Fellowship grant HST-HF2-51365.001-A awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc. for NASA under contract NAS 5-26555.
Work by CDD was performed in part under contract with the Jet Propulsion Laboratory (JPL) funded by NASA through the Sagan Fellowship Program executed by the NASA Exoplanet Science Institute.
This research used the computing resources of NERSC, a DOE Office of Science User Facility supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231.
Finally, the authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of Maunakea has always had within the indigenous Hawaiian community.  We are most fortunate to have the opportunity to conduct observations from this mountain.